home *** CD-ROM | disk | FTP | other *** search
/ The Fatted Calf / The Fatted Calf.iso / Applications / Graphics / NXPlot3d / Source / Expression.m < prev    next >
Text File  |  1994-02-19  |  28KB  |  901 lines
  1. /*
  2.     Expression.m
  3.  
  4.     The Expression class is implemented using a grammar generated by the Unix
  5.     program yacc, and a lexical scanner generated by the Unix program lex.
  6.     The only interface between these parsing modules and this class is the
  7.     function _EXPParseExpression, which actually performs the parse.  The
  8.     results of the parse are returned in two data structures: a parse tree
  9.     representing the expression and a hashtable which holds that variables
  10.     found in the expression.
  11.  
  12.     The parse tree's nodes are Term structs (declared in exprDefs.h).  There
  13.     is a node for each piece of the expression.  For example, when the
  14.     expression "A+5" is parsed, three nodes result.  The top node of the tree
  15.     represents the binary operator "+".  This top node has two subnodes.  The
  16.     first subnode represents the variable "A".  The second represents the
  17.     integer constant "5".  Since variables can appear more than once in an
  18.     expression, their nodes may occur more than once in the parse tree.
  19.     
  20.     A parse tree is evaluated by doing a depth first traversal of the tree,
  21.     bubbling up the results of each sub-tree until the final value rises
  22.     to the top.
  23.  
  24.     You may freely copy, distribute, and reuse the code in this example.
  25.     NeXT disclaims any warranty of any kind, expressed or implied, as to its
  26.     fitness for any particular use.
  27. */
  28.  
  29. #import "Graph.h"
  30.  
  31. /* declaration of methods static to this class */
  32. @interface Expression(ExpressionPrivate)
  33. - (int)_updateResults;
  34. - (EXPTermPtr)_addVarTerm:(const char *)name;
  35. @end
  36.  
  37. @implementation Expression
  38.  
  39. /* function which can be applied to parse tree nodes with applyToTerms() */
  40. typedef void ParseTreeFunc(void *info, Term *term);
  41.  
  42. static void freeContents(Expression *self);
  43. static void applyToTerms(Term *t, ParseTreeFunc *func, void *data, int mask);
  44. static void updateTerm(void *data, Term *t);
  45. static float evalTerm(Term *t, int *indices);
  46. static Term *termOfVar(Expression *self, const char *varName);
  47. static NXHashTable *makeBuiltInFuncTable(NXZone *zone);
  48. static NXHashTable *getBuiltInFuncs(void);
  49. static Function *addFuncTerm(NXHashTable *table, const char *name, int min, int max, EXPTermEvalFunc *func);
  50. static void safeFree(void **data);
  51. static EXPTermEvalFunc sinStub, cosStub, tanStub;
  52. static EXPTermEvalFunc asinStub, acosStub, atanStub;
  53. static EXPTermEvalFunc expStub, lnStub, sqrtStub, sumStub;
  54. static EXPTermEvalFunc powStub,log10Stub,jnStub,ynStub;
  55. static EXPTermEvalFunc absStub,floorStub,ceilStub,YlmStub;
  56. static char *termName(const Term *t);
  57. static void FunctionFree(const void *info, void *data);
  58. static unsigned    VarTermHash(const void *info, const void *data);
  59. static int VarTermCompare(const void *info, const void *data1, const void *data2);
  60.  
  61. /*
  62.  * Shared table of built in functions.  All Expressions which haven't had
  63.  * any application functions added to them shared this table, which contains
  64.  * just the built in functions.
  65.  */
  66. static NXHashTable *BuiltInFuncTable = NULL;
  67.  
  68. /* data for a built in function */
  69. typedef struct _BuiltInFunc {
  70.     const char *name;        /* name of the function */
  71.     EXPTermEvalFunc *func;    /* proc to call for evaluation */
  72.     int minArgs;
  73.     int maxArgs;
  74. } BuiltInFunc;
  75.  
  76. /* table of built in functions */
  77. static const BuiltInFunc FuncList[] = {
  78.     {"sin", &sinStub, 1, 1},
  79.     {"cos", &cosStub, 1, 1},
  80.     {"tan", &tanStub, 1, 1},
  81.     {"asin", &asinStub, 1, 1},
  82.     {"acos", &acosStub, 1, 1},
  83.     {"atan", &atanStub, 1, 1},
  84.     {"exp", &expStub, 1, 1},
  85.     {"ln", &lnStub, 1, 1},
  86.     {"log", &lnStub, 1, 1},
  87.     {"sqrt", &sqrtStub, 1, 1},
  88.     {"sum", &sumStub, 1, -1},
  89.     {"pow", &powStub, 2, 2},
  90.     {"log10",&log10Stub,1,1},
  91.     {"jn",&jnStub,2,2},
  92.     {"yn",&ynStub,2,2},
  93.     {"ylm",&YlmStub,4,4}, 
  94.     {"floor",&floorStub,1,1},
  95.     {"ceil",&ceilStub,1,1},
  96.     {"abs",&absStub,1,1}
  97. };
  98.  
  99. #define NUM_BUILTIN_FUNCS    (sizeof(FuncList)/sizeof(BuiltInFunc))
  100.  
  101. /* prototype used to create hashtables of variable terms */
  102. static NXHashTablePrototype VarTermProto = {&VarTermHash, &VarTermCompare, (void (*)(const void *info, void *data))&_EXPFreeTerm, 0};
  103.  
  104. - init {
  105.     [super init];
  106.     resolution = 1;
  107.     dimensions = 1;
  108.     varTerms = NXCreateHashTableFromZone(VarTermProto, 0, NULL, [self zone]);
  109.     return self;
  110. }
  111.  
  112. - free {
  113.     freeContents(self);
  114.   /* if we have a function table and its not the shared one */
  115.     if (validFuncs && validFuncs != BuiltInFuncTable)
  116.     NXFreeHashTable(validFuncs);
  117.     return [super free];
  118. }
  119.  
  120. - (BOOL)parse:(const char *)expressionString {
  121.     NXZone *zone;
  122.  
  123.     zone = [self zone];
  124.   /* clear away any results from a previous parse */
  125.     freeContents(self);
  126.     varTerms = NXCreateHashTableFromZone(VarTermProto, 0, NULL, zone);
  127.     resultsValid = NO;
  128.     text = NXCopyStringBufferFromZone(expressionString, zone);
  129.     if (!validFuncs)
  130.     validFuncs = getBuiltInFuncs();
  131.     return _EXPParseExpression(text, validFuncs, &parseTree, varTerms, zone);
  132. }
  133.  
  134. - (const char *)text {
  135.     return text;
  136. }
  137.  
  138. - setResolution:(int)count {
  139.     if (resolution != count) {
  140.       /*
  141.        * Changing the resolution means we have to recalculate next time
  142.        * we're asked for results
  143.        */
  144.     resultsValid = NO;
  145.     resolution = count;
  146.     }
  147.     return self;
  148. }
  149.  
  150. - (int)resolution {
  151.     return resolution;
  152. }
  153.  
  154. - setVar:(const char *)varName value:(float)val {
  155.     Term *t;
  156.  
  157.     t = termOfVar(self, varName);
  158.     if (!t)
  159.     t = [self _addVarTerm:varName];
  160.   /*
  161.    * If the term was previously a vector term, we change it to a variable
  162.    * term (one that has a single value).
  163.    */
  164.     if (t->tag == vectorTerm) {
  165.     t->tag = varTerm;
  166.     safeFree((void **)&t->data.vector.vals);
  167.     t->data.var.name = t->data.vector.name;
  168.     resultsValid = NO;
  169.     }
  170.     NX_ASSERT(t->tag == varTerm, "Invalid term type in setVar:value:");
  171.     if (t->data.var.val != val) {
  172.     t->data.var.val = val;
  173.     resultsValid = NO;        /* must recalc after a var is set */
  174.     }
  175.     return self;
  176. }
  177.  
  178. - (float)varValue:(const char *)varName {
  179.     Term *t;
  180.  
  181.     t = termOfVar(self, varName);
  182.     if (t) {
  183.     if (t->tag == varTerm)
  184.         return t->data.var.val;
  185.     else
  186.         NX_RAISE(expErrInvalidVarType, self, (void *)varName);
  187.     } else
  188.     NX_RAISE(expErrInvalidVarName, self, (void *)varName);
  189. }
  190.  
  191. - setVar:(const char *)varName vector:(float *)vals numVals:(int)count {
  192.     Term *t;
  193.  
  194.     resultsValid = NO;
  195.     t = termOfVar(self, varName);
  196.     if (!t)
  197.     t = [self _addVarTerm:varName];
  198.   /*
  199.    * If the term was previously a non-vector variable, we change it to a
  200.    * vector variable term.
  201.    */
  202.     if (t->tag == varTerm) {
  203.     t->tag = vectorTerm;
  204.     t->data.vector.name = t->data.var.name;
  205.     t->data.vector.hasRange = NO;
  206.     t->data.vector.vals = NULL;
  207.     t->data.vector.resolution = 0;
  208.     t->data.vector.dimension = 0;
  209.     }
  210.     NX_ASSERT(t->tag == vectorTerm, "Invalid term type in setVar:vector:");
  211.     safeFree((void **)&t->data.vector.vals);
  212.     t->data.vector.resolution = count;
  213.     t->data.vector.vals = vals;
  214.     t->data.vector.changed = YES;
  215.     resolution = count;
  216.     return self;
  217. }
  218.  
  219. - varVector:(const char *)varName vector:(float **)vals numVals:(int *)count {
  220.     Term *t;
  221.  
  222.     t = termOfVar(self, varName);
  223.     if (t) {
  224.     if (t->tag == vectorTerm) {
  225.         updateTerm(self, t);
  226.         *count = t->data.vector.resolution;
  227.         *vals = t->data.vector.vals;
  228.     } else
  229.         NX_RAISE(expErrInvalidVarType, self, (void *)varName);
  230.     } else
  231.         NX_RAISE(expErrInvalidVarName, self, (void *)varName);
  232.     return self;
  233. }
  234.  
  235. - setVar:(const char *)varName min:(float)minVal max:(float)maxVal {
  236.     Term *t;
  237.  
  238.     if (minVal > maxVal)
  239.     NX_RAISE(expErrMinMax, self, NULL);
  240.     t = termOfVar(self, varName);
  241.     if (!t)
  242.     t = [self _addVarTerm:varName];
  243.   /*
  244.    * If the term was previously a non-vector variable, we change it to a
  245.    * vector variable term.
  246.    */
  247.     if (t->tag == varTerm) {
  248.     t->tag = vectorTerm;
  249.     t->data.vector.name = t->data.var.name;
  250.     t->data.vector.hasRange = YES;
  251.     t->data.vector.vals = NULL;
  252.     t->data.vector.changed = YES;
  253.     t->data.vector.dimension = 0;
  254.     }
  255.     NX_ASSERT(t->tag == vectorTerm, "Invalid term type in setVar:min:max:");
  256.  
  257.   /*
  258.    * We optimize and do nothing if the passed values are the same as our
  259.    * current values.
  260.    */
  261.     if (t->data.vector.changed || t->data.vector.min != minVal || t->data.vector.max != maxVal) {
  262.     safeFree((void **)&t->data.vector.vals);
  263.     t->data.vector.resolution = 0;
  264.     t->data.vector.min = minVal;
  265.     t->data.vector.max = maxVal;
  266.     t->data.vector.changed = YES;
  267.     resultsValid = NO;
  268.     }    
  269.     return self;
  270. }
  271.  
  272. - var:(const char *)varName min:(float *)minVal max:(float *)maxVal {
  273.     Term *t;
  274.  
  275.     t = termOfVar(self, varName);
  276.     if (t) {
  277.     if (t->tag == vectorTerm) {
  278.         updateTerm(self, t);
  279.         *minVal = t->data.vector.min;
  280.         *maxVal = t->data.vector.max;
  281.     } else
  282.         NX_RAISE(expErrInvalidVarType, self, (void *)varName);
  283.     } else
  284.         NX_RAISE(expErrInvalidVarName, self, (void *)varName);
  285.     return self;
  286. }
  287.  
  288. - setVar:(const char *)varName dimension:(short)dimensionNum {
  289.     Term *t;
  290.  
  291.     if (dimensionNum < 0 || dimensionNum >= dimensions)
  292.     NX_RAISE(expInvalidDimension, self, (void *)varName);
  293.     resultsValid = NO;
  294.     t = termOfVar(self, varName);
  295.     if (!t)
  296.     t = [self _addVarTerm:varName];
  297.   /*
  298.    * If the term was previously a non-vector variable, we change it to a
  299.    * vector variable term.
  300.    */
  301.     if (t->tag == varTerm) {
  302.     t->tag = vectorTerm;
  303.     t->data.vector.name = t->data.var.name;
  304.     t->data.vector.hasRange = YES;
  305.     t->data.vector.vals = NULL;
  306.     t->data.vector.changed = YES;
  307.     t->data.vector.dimension = 0;
  308.     }
  309.     NX_ASSERT(t->tag == vectorTerm, "Invalid term type in setVar:dimension:");
  310.     t->data.vector.dimension = dimensionNum;
  311.     return self;
  312. }
  313.  
  314. - var:(const char *)varName dimension:(short *)dimensionNum {
  315.     Term *t;
  316.  
  317.     t = termOfVar(self, varName);
  318.     if (t) {
  319.     if (t->tag == vectorTerm)
  320.         *dimensionNum = t->data.vector.dimension;
  321.     else
  322.         NX_RAISE(expErrInvalidVarType, self, (void *)varName);
  323.     } else
  324.         NX_RAISE(expErrInvalidVarName, self, (void *)varName);
  325.     return self;
  326. }
  327.  
  328. - (float)resultValue {
  329.     [self _updateResults];
  330.     return *results;
  331. }
  332.  
  333. - resultsVector:(float **)vals numVals:(int *)count {
  334.     *count = [self _updateResults];
  335.     *vals = results;
  336.     return self;
  337. }
  338.  
  339. - resultsMin:(float *)minVal max:(float *)maxVal {
  340.     [self _updateResults];
  341.     *minVal = resultsMin;
  342.     *maxVal = resultsMax;
  343.     return self;
  344. }
  345.  
  346. - setDimensions:(short)count {
  347.     if (dimensions != count)
  348.     resultsValid = NO;
  349.     dimensions = count;
  350.     return self;
  351. }
  352.  
  353. - (short)dimensions {
  354.     return dimensions;
  355. }
  356.  
  357. /*
  358.  * Since the varables are all in a NXHashTable, we just use the NXHashTable
  359.  * functions to enumerate through the names of the variables.
  360.  */
  361.  
  362. - (EXPEnumState)beginVariableEnumeration {
  363.     NXHashState *state;
  364.  
  365.     state = NXZoneMalloc([self zone], sizeof(NXHashState));
  366.     *state = NXInitHashState(varTerms);
  367.     return state;
  368. }
  369.  
  370. - (const char *)nextVariable:(EXPEnumState)state {
  371.     const Term *t;
  372.  
  373.     if (NXNextHashState(varTerms, (NXHashState *)state, (void **)&t))
  374.     return termName(t);
  375.     else
  376.     return NULL;
  377. }
  378.  
  379. - (void)endVariableEnumeration:(EXPEnumState)state {
  380.     free(state);
  381. }
  382.  
  383. - addFuncTerm:(const char *)name minArgs:(int)min maxArgs:(int)max
  384.                     evalFunc:(EXPTermEvalFunc *)func {
  385.     Function *existingType;
  386.  
  387.   /*
  388.    * If we dont have a function table yet, get one.  If we do have one, but
  389.    * its the shared built in table, get a new one that we can safely modify.
  390.    */
  391.     if (!validFuncs || validFuncs == BuiltInFuncTable)
  392.     validFuncs = makeBuiltInFuncTable([self zone]);
  393.     existingType = addFuncTerm(validFuncs, name, min, max, func);
  394.     if (existingType)
  395.     NX_RAISE(expFuncTypeInUse, self, existingType);
  396.     return self;
  397. }
  398.  
  399. - removeFuncTerm:(const char *)name {
  400.     Function *realFunc;
  401.     Function key;
  402.  
  403.     if (validFuncs) {
  404.       /* look up the func term by name in our table of functions */
  405.     key.name = (char *)name;
  406.     if (realFunc = NXHashGet(validFuncs, &key)) {
  407.       /*
  408.        * If we are using the shared table of built ins, get a new table
  409.        * that we can safely modify (this covers the case of someone
  410.        * wishing to remove a built in function, possibly to redefine it).
  411.        */
  412.         if (validFuncs == BuiltInFuncTable) {
  413.         validFuncs = makeBuiltInFuncTable([self zone]);
  414.         realFunc = NXHashGet(validFuncs, &key);
  415.         }
  416.         NXHashRemove(validFuncs, realFunc);
  417.         free(realFunc);
  418.         return self;
  419.     }
  420.     }
  421.     return self;
  422. }
  423.  
  424. /*
  425.  * This is a utililty routine that recursively applies a function to all
  426.  * terms of a parse tree. data is a blind pointer that is simply passed
  427.  * along through the function call.  The function is only called on terms
  428.  * whose type match the given mask.
  429.  */
  430. static void applyToTerms(Term *t, ParseTreeFunc *func, void *data, int mask) {
  431.     int i;
  432.     Term **tPtr;
  433.  
  434.     for (i = t->numSubterms, tPtr = t->subterms; i--; tPtr++)
  435.     applyToTerms(*tPtr, func, data, mask);
  436.     if (t->tag & mask)
  437.     (*func)(data, t);
  438. }
  439.  
  440. /*
  441.  * Empties the contents of the Expression.  Since the variable terms can
  442.  * exist multiple times in the tree, we first run through the tree and
  443.  * free all the nodes except that variable nodes.  Then we free the hash
  444.  * table of variable terms, including the terms themselves. 
  445.  */
  446. static void freeContents(Expression *self) {
  447.     safeFree((void **)&self->text);
  448.   /* free the non-variable terms */
  449.     if (self->parseTree)
  450.     applyToTerms(self->parseTree, _EXPFreeTerm, NULL, ~(varTerm|vectorTerm));
  451.   /* free the shared variable terms */
  452.     NXFreeHashTable(self->varTerms);
  453.     safeFree((void **)&self->results);
  454. }
  455.  
  456. /*
  457.  * Allocates a new term of the given type, with room for subterms.  The
  458.  * subterms themselves follow as a variable number of arguments.  The are
  459.  * copied into the subterms list of the new term.
  460.  */
  461. Term *_EXPAllocTerm(NXZone *zone, TermTag tag, int numSubterms, ...) {
  462.     Term *t;
  463.     int i;
  464.     va_list args;
  465.  
  466.     t = NXZoneCalloc(zone, sizeof(Term) + (numSubterms-1) * sizeof(Term *), 1);
  467.     t->tag = tag;
  468.     t->numSubterms = numSubterms;
  469.     va_start(args, numSubterms);
  470.     for (i = 0; i < numSubterms; i++)
  471.     t->subterms[i] = va_arg(args, Term *);
  472.     va_end(args);
  473.     return t;
  474. }
  475.  
  476. /*
  477.  * Frees a term and any associated data.  This routine can be used as the
  478.  * free function of a NXHashTable prototype, or as a proc passed to
  479.  * applyToTerms().
  480.  */
  481. void _EXPFreeTerm(void *info, Term *data) {
  482.     Term *t = (Term *)data;
  483.  
  484.     if (t->tag == vectorTerm) {
  485.     free(t->data.vector.vals);
  486.     free(t->data.vector.name);
  487.     } else if (t->tag == varTerm)
  488.     free(t->data.var.name);
  489.     free(t);
  490. }
  491.  
  492. /*
  493.  * Makes sure a term is up to date.  Since terms recalculate any internal
  494.  * state lazily, this must be called before making use of a term's value.
  495.  * We apply this function recursively to all terms before evaluating an
  496.  * Expression, and apply to any term if we are asked to return the values
  497.  * held within the term (for example, from -varVector:vector:numVals:).
  498.  */
  499. static void updateTerm(void *data, Term *t) {
  500.     Expression *self = data;
  501.  
  502.     if (t->tag != vectorTerm)
  503.     return;        /* only vector terms require work to stay up to date */
  504.  
  505.   /* Ensure this term has the same resolution as the rest of the Expression. */ 
  506.     if (self->resolution != t->data.vector.resolution) {
  507.       /*
  508.        * We can change its resolution if we interpolate values for this term
  509.        * within a range.  Else if we were passed a list of values for this
  510.        * term, its an exception if the resolution is no longer in sync.
  511.        */
  512.     if (t->data.vector.hasRange) {
  513.         safeFree((void **)&t->data.vector.vals);
  514.         t->data.vector.vals = NXZoneMalloc([self zone],
  515.                     self->resolution * sizeof(float));
  516.         t->data.vector.resolution = self->resolution;
  517.         t->data.vector.changed = YES;    /* remember to re-interpolate */
  518.     } else
  519.         NX_RAISE(expErrResolutionMismatch, self,
  520.                     (void *)t->data.vector.name);
  521.     }
  522.  
  523.     if (t->data.vector.changed) {
  524.     if (t->data.vector.hasRange) {
  525.       /* interpolate a list of values between min and max */
  526.         int i;
  527.         float delta;
  528.         float *val, *prevVal;
  529.  
  530.         i = self->resolution - 1;
  531.         if (i) {
  532.         delta = (t->data.vector.max - t->data.vector.min) / i;
  533.         prevVal = t->data.vector.vals;
  534.         *prevVal = t->data.vector.min;
  535.         val = prevVal + 1;
  536.         while (i--)
  537.             *val++ = *prevVal++ + delta;
  538.         *(val-1) = t->data.vector.max;    /* to be sure we hit max */
  539.         } else
  540.         *t->data.vector.vals = t->data.vector.min;
  541.     } else {
  542.       /* scan the list of values passed in to find the min and max */
  543.         int i;
  544.         float newMin, newMax;
  545.         float *val;
  546.  
  547.         val = t->data.vector.vals;
  548.         newMin = newMax = *val++;
  549.         for (i = t->data.vector.resolution; i--; val++) {
  550.         if (*val > newMax)
  551.             newMax = *val;
  552.         else if (*val < newMin)
  553.             newMin = *val;
  554.         }
  555.         t->data.vector.min = newMin;
  556.         t->data.vector.max = newMax;
  557.     }
  558.     t->data.vector.changed = NO;    /* we're now up to date */
  559.     }
  560. }
  561.  
  562. /*
  563.  * Ensures that the results calculated for this Expression are up to date.
  564.  * We first make sure all the terms in the parse tree are up to date by
  565.  * applying updateTerm() to all of them.  We then evaluate the Expression
  566.  * for every n times, depending on the resolution, storing all the results
  567.  * and calculating their min and max. 
  568.  */
  569. - (int)_updateResults {
  570. #define MAX_INDICES 10
  571.     int i, j;
  572.     float *f;
  573.     int indicesBuffer[MAX_INDICES];
  574.     int *indices = indicesBuffer;
  575.     int totalVals;
  576.     NXZone *zone = [self zone];
  577.  
  578.     totalVals = resolution;
  579.     for (i = 1; i < dimensions; i++)
  580.     totalVals *= resolution;
  581.     if (!resultsValid) {
  582.     if (parseTree) {
  583.         applyToTerms(self->parseTree, updateTerm, self, -1);
  584.         safeFree((void **)&results);
  585.         results = NXZoneMalloc(zone, sizeof(float) * totalVals);
  586.         if (dimensions > MAX_INDICES)
  587.         indices = NXZoneMalloc(zone, sizeof(int) * dimensions);
  588.         bzero(indices, sizeof(int) * dimensions);
  589.         *results = evalTerm(parseTree, indices);
  590.         resultsMin = *results;
  591.         resultsMax = *results;
  592.         for (i = 1, f = results + 1; i < totalVals; i++, f++) {
  593.  
  594.           /*
  595.            * Increment the indices.  The first one always changes.  The
  596.            * others change only when the previous one wraps around (like
  597.            * an odometer with base = resolution).
  598.            */
  599.         for (j = 0; j < dimensions; j++) {
  600.             indices[j] = (indices[j] + 1) % resolution;
  601.             if (indices[j])
  602.             break;
  603.         }
  604.  
  605.           /*
  606.            * We pass the loop indices down through the evaluation recursion
  607.            * so vector terms can know which element of their vectors they
  608.            * should use for this evaluation.
  609.            */
  610.         *f = evalTerm(parseTree, indices);
  611.         if (*f > resultsMax)
  612.             resultsMax = *f;
  613.         else if (*f < resultsMin)
  614.             resultsMin = *f;
  615.         }
  616.         if (dimensions > MAX_INDICES)
  617.         NXZoneFree(zone, indices);
  618.     } else
  619.         NX_RAISE(expErrNoText, self, NULL);
  620.     resultsValid = YES;
  621.     }
  622.     return totalVals;
  623. }
  624.  
  625. #define MAX_ARGS    50
  626.  
  627. /*
  628.  * Evaluates a particular term.  In order to evaluate itself, any term with
  629.  * subterms must recursively evaluate those first.
  630.  */
  631. static float evalTerm(Term *t, int *indices) {
  632.     float result = 0.0;        /* initted to quiet -Wall */
  633.     float base;
  634.     float argsBuffer[MAX_ARGS];
  635.     float *args;
  636.     int i;
  637.  
  638.     switch (t->tag) {
  639.     case constantTerm:
  640.         NX_ASSERT(t->numSubterms == 0, "Wrong #subterms in evalTerm");
  641.         result = t->data.constant.val;
  642.         break;
  643.     case varTerm:
  644.         NX_ASSERT(t->numSubterms == 0, "Wrong #subterms in evalTerm");
  645.         result = t->data.var.val;
  646.         break;
  647.     case vectorTerm:
  648.         NX_ASSERT(t->numSubterms == 0, "Wrong #subterms in evalTerm");
  649.         result = t->data.vector.vals[indices[t->data.vector.dimension]];
  650.         break;
  651.     case binOpTerm:
  652.         NX_ASSERT(t->numSubterms == 2 ||
  653.             (t->data.binOp.op == '-' && t->numSubterms == 1),
  654.                     "Wrong #subterms in evalTerm");
  655.         switch (t->data.binOp.op) {
  656.         case '+':
  657.             result = evalTerm(t->subterms[0], indices) +
  658.                 evalTerm(t->subterms[1], indices);
  659.             break;
  660.         case '-':
  661.             if (t->numSubterms == 2)
  662.             result = evalTerm(t->subterms[0], indices) -
  663.                     evalTerm(t->subterms[1], indices);
  664.             else
  665.             result = - evalTerm(t->subterms[0], indices);
  666.             break;
  667.         case '*':
  668.             result = evalTerm(t->subterms[0], indices) *
  669.                 evalTerm(t->subterms[1], indices);
  670.             break;
  671.         case '/':
  672.             result = evalTerm(t->subterms[0], indices) /
  673.                 evalTerm(t->subterms[1], indices);
  674.             break;
  675.         case '%':
  676.             result = (int)rint(evalTerm(t->subterms[0], indices)) %
  677.                 (int)rint(evalTerm(t->subterms[1], indices));
  678.             break;
  679.         case '^':
  680.           /* optimize for raising to an integral power */
  681.             if (t->subterms[1]->tag == constantTerm &&
  682.                 t->subterms[1]->data.constant.isInt &&
  683.                 t->subterms[1]->data.constant.val >= 1) {
  684.             result = base = evalTerm(t->subterms[0], indices);
  685.             for (i = t->subterms[1]->data.constant.val; --i; )
  686.                 result *= base;
  687.             } else
  688.             result = pow(evalTerm(t->subterms[0], indices),
  689.                     evalTerm(t->subterms[1], indices));
  690.             break;
  691.         default:
  692.             NX_ASSERT(FALSE, "Unknown binary op type in evalTerm");
  693.         }
  694.         break;
  695.     case funcTerm:
  696.       /*
  697.        * For functions, we first ensure we have a large enough buffer
  698.        * for the values of all the arguments.  If there are few enough
  699.        * arguments, we use a buffer on the stack instead of thrashing
  700.        * the heap.  We then buffer up all the results of evaluating
  701.        * the arguments, and then pass this array of argument values
  702.        * to the proc that we use to evaluate this type of function.
  703.        */
  704.         if (t->numSubterms > MAX_ARGS)
  705.         args = NXZoneMalloc(NXDefaultMallocZone(),
  706.                     sizeof(float) * t->numSubterms);
  707.         else
  708.         args = argsBuffer;
  709.         for (i = 0; i < t->numSubterms; i++)
  710.         args[i] = evalTerm(t->subterms[i], indices);
  711.         result = t->data.func.type->evalFunc(t->numSubterms, args);
  712.             if (t->numSubterms > MAX_ARGS)
  713.         NXZoneFree(NXDefaultMallocZone(), args);
  714.         break;
  715.     default:
  716.         NX_ASSERT(FALSE, "Invalid term type in evalTerm");
  717.     }
  718.     return result;
  719. }
  720.  
  721. /* Utility routine to look up a variable by name in the variable hashtable. */
  722. static Term *termOfVar(Expression *self, const char *varName) {
  723.     Term key;
  724.  
  725.     if (self->parseTree) {
  726.     key.tag = varTerm;
  727.     key.data.var.name = (char *)varName;
  728.     return NXHashGet(self->varTerms, &key);
  729.     } else
  730.     NX_RAISE(expErrNoText, self, NULL);
  731. }
  732.  
  733. /* adds a variable term to the Expressions hashtable of them */
  734. - (Term *)_addVarTerm:(const char *)name {
  735.     Term *newTerm;
  736.     Term *existingTerm;
  737.  
  738.     newTerm = _EXPAllocTerm([self zone], varTerm, 0);
  739.     newTerm->data.var.name = NXCopyStringBufferFromZone(name, [self zone]);
  740.     existingTerm = NXHashInsertIfAbsent(varTerms, newTerm);
  741.     NX_ASSERT(existingTerm == newTerm, "_addVarTerm: called with existing term");
  742.     return newTerm;
  743. }
  744.  
  745. /* frees some storage, NULL'ing out the pointer */
  746. static void safeFree(void **data) {
  747.     free(*data);
  748.     *data = NULL;
  749. }
  750.  
  751. /* free function used in the NXHashTable prototype for functions */
  752. static void FunctionFree(const void *info, void *data) {
  753.     free(((Function *)data)->name);
  754.     free(data);
  755. }
  756.  
  757. /*
  758.  * Adds a func term to a HashTable of them.  Returns any existing entry
  759.  * with the same name, else NULL.
  760.  */
  761. static Function *addFuncTerm(NXHashTable *table, const char *name, int min, int max, EXPTermEvalFunc *func) {
  762.     Function *newFunc;
  763.     Function *existingType;
  764.  
  765.     newFunc = NXZoneMalloc(NXZoneFromPtr(table), sizeof(Function));
  766.     newFunc->name = NXCopyStringBufferFromZone(name, NXZoneFromPtr(table));
  767.     newFunc->minArgs = min;
  768.     newFunc->maxArgs = max;
  769.     newFunc->evalFunc = func;
  770.     existingType = NXHashInsertIfAbsent(table, newFunc);
  771.     if (existingType != newFunc)
  772.     return existingType;
  773.     else
  774.     return NULL;
  775. }
  776.  
  777. /*
  778.  * Returns a global table of all built in functions.  This table is shared
  779.  * by expressions that dont have application functions added to them.
  780.  */
  781. static NXHashTable *getBuiltInFuncs(void) {
  782.     if (!BuiltInFuncTable)
  783.     BuiltInFuncTable = makeBuiltInFuncTable(NXDefaultMallocZone());
  784.     return BuiltInFuncTable;
  785. }
  786.  
  787. /* Returns a new hashtable of all built in functions. */
  788. static NXHashTable *makeBuiltInFuncTable(NXZone *zone) {
  789.     NXHashTable *table;
  790.     NXHashTablePrototype FuncTermProto;
  791.     const BuiltInFunc *bif;
  792.     int i;
  793.  
  794.     FuncTermProto = NXStrStructKeyPrototype;
  795.     FuncTermProto.free = FunctionFree;
  796.     table = NXCreateHashTableFromZone(FuncTermProto, NUM_BUILTIN_FUNCS,
  797.                                 NULL, zone);
  798.     for (i = NUM_BUILTIN_FUNCS, bif = FuncList; i--; bif++)
  799.     (void)addFuncTerm(table, bif->name, bif->minArgs, bif->maxArgs, bif->func);
  800.     return table;
  801. }
  802.  
  803. /* Returns the name of a term. */
  804. static char *termName(const Term *t) {
  805.     switch (t->tag) {
  806.     case varTerm:
  807.         return t->data.var.name;
  808.     case vectorTerm:
  809.         return t->data.vector.name;
  810.     default:
  811.         NX_ASSERT(FALSE, "Bogus term type in VarTermHash");
  812.         return NULL;
  813.     }
  814. }
  815.  
  816. /* hashing function for variable terms.  Used in hashtable prototypes. */
  817. static unsigned    VarTermHash(const void *info, const void *data) {
  818.     return NXStrHash(info, termName(data));
  819. }
  820.  
  821. /* comparison function for variable terms.  Used in hashtable prototypes. */
  822. static int VarTermCompare(const void *info, const void *data1, const void *data2) {
  823.     return NXStrIsEqual(info, termName(data1), termName(data2));
  824. }
  825.  
  826. @end
  827.  
  828. /* These procs implement the built in functions */
  829.  
  830. static float sinStub(int numArgs, float *arg)    { return sin(*arg); }
  831. static float cosStub(int numArgs, float *arg)    { return cos(*arg); }
  832. static float tanStub(int numArgs, float *arg)    { return tan(*arg); }
  833. static float asinStub(int numArgs, float *arg)    { return asin(*arg); }
  834. static float acosStub(int numArgs, float *arg)    { return acos(*arg); }
  835. static float atanStub(int numArgs, float *arg)    { return atan(*arg); }
  836. static float expStub(int numArgs, float *arg)    { return exp(*arg); }
  837. static float lnStub(int numArgs, float *arg)    { return log(*arg); }
  838. static float sqrtStub(int numArgs, float *arg)    { return sqrt(*arg); }
  839.  
  840. static float sumStub(int numArgs, float *arg) {
  841.     float sum = 0.0;
  842.  
  843.     while (numArgs--)
  844.     sum += *arg++;
  845.     return sum;
  846. }
  847.  
  848. /* added 2/94  Steve Ludtke */
  849. static float powStub(int numArgs, float *arg)     { return pow(arg[0],arg[1]); }
  850. static float log10Stub(int numArgs, float *arg) { return log10(arg[0]); }
  851. static float jnStub(int numArgs, float *arg)     
  852.     { return jn((int)arg[0],arg[1]); }
  853. static float ynStub(int numArgs, float *arg)     
  854.     { return yn((int)arg[0],arg[1]); }
  855. static float absStub(int numArgs, float *arg)    { return fabs(arg[0]); }
  856. static float floorStub(int numArgs, float *arg)    { return floor(arg[0]); }
  857. static float ceilStub(int numArgs, float *arg)    { return ceil(arg[0]); }
  858.  
  859.  
  860. static float YlmStub(int numArgs, float *arg) 
  861. {
  862. int l,m;
  863. float theta,phi;
  864.  
  865. if (numArgs!=4) return(0.0);
  866. l=(int)arg[0];
  867. m=(int)arg[1];
  868. theta=arg[2];
  869. phi=arg[3];
  870.  
  871. switch(l) {
  872.     case 0:
  873.         if (m==0) return(sqrt(1.0/(4.0*M_PI)));
  874.         break;
  875.     case 1:
  876.         if (m==-1) return(1.0/2.0*sqrt(1.5/M_PI)*sin(theta)*sin(phi));
  877.         if (m==0) return(1.0/2.0*sqrt(3.0/M_PI)*cos(theta));
  878.         if (m==1) return(-1.0/2.0*sqrt(1.5/M_PI)*sin(theta)*cos(phi));
  879.         break;
  880.     case 2:
  881.         if (m==-2) return(1.0/4.0*sqrt(15.0/(2.0*M_PI))*pow(sin(theta),2.0)*sin(2.0*phi));
  882.         if (m==-1) return(1.0/2.0*sqrt(15.0/(2.0*M_PI))*sin(theta)*cos(theta)*sin(phi));
  883.         if (m==0) return(1/4.0*sqrt(5.0/M_PI)*(3.0*pow(cos(theta),2.0)-1.0));
  884.         if (m==1) return(-1.0/2.0*sqrt(15.0/(2.0*M_PI))*sin(theta)*cos(theta)*cos(phi));
  885.         if (m==2) return(1.0/4.0*sqrt(15.0/(2.0*M_PI))*pow(sin(theta),2.0)*cos(2.0*phi));
  886.         break;
  887.     case 3:
  888.         if (m==-3) return(1.0/8.0*sqrt(35.0/M_PI)*pow(sin(theta),3.0)*sin(3.0*phi));
  889.         if (m==-2) return(1.0/4.0*sqrt(52.5/M_PI)*pow(sin(theta),2.0)*cos(theta)*sin(2.0*phi));
  890.         if (m==-1) return(1.0/8.0*sqrt(21.0/M_PI)*sin(theta)*(5.0*pow(cos(theta),2.0)-1.0)*sin(phi));
  891.         if (m==0) return(1.0/4.0*sqrt(7/M_PI)*(5.0*pow(cos(theta),3.0)-3.0*cos(theta)));
  892.         if (m==1) return(-1.0/8.0*sqrt(21.0/M_PI)*sin(theta)*(5.0*pow(cos(theta),2.0)-1.0)*cos(phi));
  893.         if (m==2) return(1.0/4.0*sqrt(52.5/M_PI)*pow(sin(theta),2.0)*cos(theta)*cos(2.0*phi));
  894.         if (m==3) 
  895. return(-1.0/8.0*sqrt(35.0/M_PI)*pow(sin(theta),3.0)*cos(3.0*phi));
  896.         break;
  897. }
  898.  
  899. return(0.0);
  900. }
  901.